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Abstract 

We  present  a  number  of  complete  sets  of  series  expansion  formulas  for  neutrino  oscillation 
probabilities  in  matter  of  constant  density  for  three  flavors.  In  particular,  we  study  expansions 
in  the  mass  hierarchy  parameter  a  =  and  mixing  parameter  S13  =  sin^ia  up  to 

second  order  and  expansions  only  in  a  and  only  in  S13  up  to  first  order.  For  each  type  of 
expansion  we  also  present  the  corresponding  formulas  for  neutrino  oscillations  in  vacuum.  We 
perform  a  detailed  analysis  of  the  accuracy  of  the  different  sets  of  series  expansion  formulas 
and  investigate  which  type  of  expansion  is  most  accurate  in  different  regions  of  the  parameter 
space  spanned  by  the  neutrino  energy  E,  the  baseline  length  L,  and  the  expansion  parameters 
a  and  S13.  We  also  present  the  formulas  for  series  expansions  in  a  and  in  S13  up  to  first  order 
for  the  case  of  arbitrary  matter  density  profiles.  Furthermore,  it  is  shown  that  in  general  all 
the  18  neutrino  and  antineutrino  oscillation  probabilities  can  be  expressed  through  just  two 
independent  probabilities. 
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1  Introduction 


The  discovery  of  neutrino  oscillations  in  atmospheric,  solar,  and  reactor  neutrino  experiments 
has  turned  neutrino  physics  into  one  of  the  most  exciting  and  active  helds  of  particle  physics. 
By  now,  a  signihcant  amount  of  information  on  neutrino  properties  has  been  obtained.  The 
results  of  the  atmospheric  neutrino  experiments  [1-3]  and  the  K2K  accelerator  neutrino  ex¬ 
periment  [4]  have  allowed  the  determination  of  the  fundamental  neutrino  parameters 
and  6*23  fo  an  accuracy  of  about  30%  and  15%,  respectively,  while  the  solar  neutrino  experi¬ 
ments  [5]  and  the  KamLAND  reactor  neutrino  experiment  [6]  have  measured  the  parameters 
and  6*12  to  an  accuracy  of  about  15%.  In  fact,  this  means  that  neutrino  physics  is 
now  entering  an  era  of  precision  measurements  of  the  neutrino  oscillation  parameters.  Future 
experiments  with  superbeams  and  neutrino  factories  will  determine  the  “atmospheric”  and 
“solar”  neutrino  oscillation  parameters  to  an  accuracy  of  the  order  of  1%.  They  are  also 
expected  to  measure  the  elusive  leptonic  mixing  angle  ^13  for  which  at  present  only  an  upper 
limit  exists  [7],  or  to  put  a  more  stringent  limit  on  this  angle,  as  well  as  to  clarify  the  issues 
of  the  neutrino  mass  hierarchy  and  possibly  of  leptonic  CP  violation.  Important  information 
on  neutrino  oscillation  parameters  can  also  be  obtained  from  the  future  atmospheric,  solar, 
reactor,  and  supernova  neutrino  experiments. 

Increasing  accuracy  and  reach  of  the  present  and  especially  of  forthcoming  experiments 
put  forward  new  and  more  challenging  demands  to  the  theoretical  description  of  neutrino 
oscillations.  In  order  to  be  able  to  determine  the  fundamental  neutrino  parameters  from  the 
data  with  high  precision,  one  needs,  among  other  things,  very  accurate  theoretical  expressions 
for  the  probabilities  of  neutrino  oscillations  in  matter  and  in  vacuum.  While  these  proba¬ 
bilities  can,  in  principle,  be  calculated  numerically  with  any  requisite  accuracy,  it  is  highly 
desirable  to  have  also  analytic  expressions  for  them.  Such  analytic  expressions  would  reveal 
the  basic  dependence  of  the  neutrino  oscillation  probabilities  on  the  fundamental  neutrino 
parameters  and  on  the  characteristics  of  the  experiment  and  thus  facilitate  the  choice  of  the 
experimental  setup  as  well  as  the  analysis  of  the  data.  They  would  also  help  to  understand 
the  physics  underlying  various  flavor  transitions  and  to  resolve  the  parameter  degeneracies 
and  other  ambiguities,  such  as  fundamental  versus  matter-induced  CP  violation. 

The  purpose  of  this  paper  is  to  present  a  collection  of  approximate  analytic  formulas 
for  the  neutrino  oscillation  probabilities.  We  derive  a  number  of  complete  sets  of  series 
expansion  formulas  for  the  three-flavor  neutrino  oscillation  probabilities  in  matter  and  in 
vacuum.  The  probabilities  are  expanded  in  the  mass  hierarchy  parameter  a  =  Arn^i/ Aril'll, 
the  mixing  parameter  S13  =  sin  6*13,  or  in  both  of  them.  We  also  study  in  detail  the  accuracy 
of  the  obtained  expressions  in  different  regions  of  the  parameter  space  and  identify  the  “best 
choice”  in  each  case  of  interest. 

Before  proceeding  to  present  our  results,  we  give  here  a  brief  overview  of  the  previous 
work  on  the  subject.  Analytic  formulas  for  three-flavor  neutrino  oscillation  probabilities  have 
been  derived  in  a  number  of  papers.  Exact  formulas  for  the  neutrino  oscillation  probabilities 
in  vacuum  can  be  found,  e.g.,  in  Ref.  [8].  Exact  expressions  can  also  be  obtained  in  the 
case  of  three-flavor  neutrino  oscillations  in  matter  of  constant  density  [9-16].  However,  the 
corresponding  formulas  are  rather  complicated  and  not  easily  tractable.  This  also  applies 
to  the  exact  analytic  three-flavor  formulas  obtained  for  some  special  cases  of  non-uniform 
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matter  density:  linear  matter  density  [17]  and  exponentially  varying  matter  density  [18]. 

Approximate  analytic  formulas  for  three-flavor  neutrino  oscillation  probabilities  in  mat¬ 
ter  have  been  derived  in  a  number  of  papers,  see  e.g.,  Refs.  [19-24].  In  Refs.  [19,20],  the 
Mikheyev-Smirnov-Wolfenstein  (MSW)  resonances  [25]  in  matter  of  varying  density  were 
studied  assuming  large  separation  and  independence  of  the  high-density  and  low-density  res¬ 
onances.  In  Refs.  [20,21,23],  the  adiabatic  approximation  was  used  to  derive  three-flavor  neu¬ 
trino  oscillation  probabilities  in  matter  of  varying  density,  whereas  in  Ref.  [24],  the  Magnus 
expansion  for  the  time  evolution  operator  was  used,  which  is  equivalent  to  the  anti-adiabatic 
approximation. 

In  a  number  of  papers  an  approach  similar  to  ours  was  adopted,  i.e.,  the  neutrino  os¬ 
cillation  probabilities  were  expanded  in  the  small  parameters  a,  S13,  or  in  both  of  them. 
In  Ref.  [26],  exact  analytic  expressions  for  three-flavor  neutrino  oscillation  probabilities  in 
matter  with  an  arbitrary  density  prohle  were  obtained  in  the  limit  a  — 0  by  reducing  the 
problem  to  an  effective  two-flavor  one.  In  Ref.  [27],  a  similar  approach  was  employed  to  ob¬ 
tain  the  neutrino  oscillation  probabilities  in  matter  of  arbitrarily  varying  density  in  the  limit 
^13  0,  whereas  in  Refs.  [28,29],  expressions  up  to  hrst  order  in  S13  were  derived.  For  the 

case  of  matter  of  constant  density,  the  limit  a  ^  0  was  considered  in  Ref.  [30].  Expansions 
up  to  hrst  order  in  a  were  carried  out  in  Refs.  [31,32].  In  Refs.  [33,34],  both  the  solar  mass 
squared  difference  Am^i  and  matter  effects  were  treated  as  perturbations  and  the  transition 
probabilities  up  to  hrst  order  in  them  were  derived.  Expansions  in  both  a  and  S13  were  used 
in  Refs.  [35-38].  We  note  that  the  neutrino  oscillation  probabilities  derived  in  most  of  the 
above-mentioned  papers  either  did  not  constitute  a  complete  set  {i.e.,  a  set  of  probabilities 
from  which  the  probabilities  in  all  channels  can  be  obtained),  or  contained  expressions  from 
which  some  terms  were  missing,  especially  in  the  case  of  the  probabilities  P{i^fj,  — ^  and 
P{^ti  ^t)-  Thus,  to  the  best  of  our  knowledge,  our  study  is  the  hrst  one  in  which,  for  the 
case  of  matter  of  constant  density,  complete  and  consistent  expansions  in  a  and  S13  to  second 
order  and  expansions  only  in  a  and  only  in  S13  up  to  hrst  order  are  performed. 

The  paper  is  organized  as  follows.  In  Sec.  2  we  set  the  general  formalism  and  notation, 
discuss  some  general  relations  satished  by  the  oscillation  probabilities,  and  show  that  all 
the  18  neutrino  and  antineutrino  oscillation  probabilities  can  be  expressed  through  just  two 
independent  probabilities.  In  Sec.  3  we  present  a  series  expansion  of  the  neutrino  oscillation 
probabilities  up  to  second  order  in  both  a  and  S13  (the  so-called  “double  expansion”),  whereas 
in  Secs.  4  and  5  we  consider  the  probabilities  expanded  up  to  hrst  order  in  a  and  S13, 
respectively  (“single  expansions”).  In  all  three  cases,  we  give  the  probabilities  for  matter  of 
constant  density  (Secs.  3.1,  4.1,  and  5.1)  and  in  vacuum  (Secs.  3.2,  4.2,  and  5.2).  We  also 
compare  our  formulas  with  the  corresponding  expressions  existing  in  the  literature  (when 
available),  pointing  out  agreements  and  disagreements.  In  Sec.  6,  we  discuss  the  qualitative 
behavior  of  the  neutrino  oscillation  probabilities,  the  relevance  of  matter  ehects,  and  give  a 
detailed  evaluation  of  the  accuracy  of  the  various  formulas.  Furthermore,  we  comment  on  the 
application  of  our  probability  formulas  to  neutrino  oscillation  experiments.  We  summarize 
our  results  in  Sec.  7.  Finally,  several  methods  that  have  been  used  to  derive  the  formulas  are 
presented  in  the  appendices.  In  Appendix  A,  we  describe  the  perturbative  diagonalization 
of  the  ehective  Hamiltonian  of  the  neutrino  system  in  the  case  of  matter  of  constant  density, 
while  in  Appendix  B,  the  details  of  the  perturbative  expansion  of  the  neutrino  evolution 
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equation  in  the  case  of  arbitrary  matter  density  profiles  are  given. 

2  General  formalism  and  notation 

We  consider  three-flavor  neutrino  oscillations  and  adopt  the  standard  parameterization  of  the 
leptonic  mixing  matrix  U  [39]: 

U  =  023Us0nUl0u 

/  C12C13  S12C13 

=  --S12C23  -  Ci2Si3S23e“^'=P  C12C23  -  Sl2-Sl3S23e‘'^'^P  C13S23  •  (1) 

\  S12S23  -  Ci2Sl3C23e''^‘^^  -C12S23  -  -Si2Sl3C23e‘'^'=P  C13C23  / 

Here  Oij  is  the  orthogonal  rotation  matrix  in  the  ij-plane  which  depends  on  the  mixing 
angle  Us  =  diag(l,  1,  e*'^'^’’),  (5cp  being  the  Dirac-type  CP-violating  phase,  Sij  =  sin  9ij 
and  Cij  =  cos  6ij.  In  the  three-flavor  case  there  are  also,  in  general,  two  Majorana-type  CP- 
violating  phases;  however,  these  phases  do  not  affect  neutrino  oscillations,  and  therefore  will 
not  be  considered  here.  Without  loss  of  generality,  one  can  assume  all  the  mixing  angles  to 
he  in  the  first  quadrant  (he.,  between  0  and  7i/2),  while  the  CP-violating  phase  (5cp  is  allowed 
to  lie  in  the  interval  [0,  27r]. 

Let  us  denote  by  Pa^  =  P{i'a  ^13)  the  probability  of  transition  from  a  neutrino  flavor 
a  to  a  neutrino  flavor  (3,  and  similarly  for  antineutrino  flavors,  he.,  P^^  =  P{^a  ^p)- 
In  general,  the  three-flavor  neutrino  oscillation  probabilities  in  matter  Pap  depend  on  eight 
parameters  and  one  function: 

Pap  =  Pap{^ml^,^ml^,0i2, 013,023, Scp;E,L,V{x)),  a,(3  =  e,n,T.  (2) 

Here  Amh  =  mj  —  are  the  neutrino  mass  squared  differences,  E  is  the  neutrino  en¬ 
ergy,  L  is  the  baseline  length,  and  V (x)  is  the  matter-induced  effective  potential,  x  E  [0,  L] 
being  the  coordinate  along  the  neutrino  path.  The  neutrino  mass  squared  differences,  the 
leptonic  mixing  angles,  and  the  CP-violating  phase  are  fundamental  parameters  and  thus 
experiment-independent,  whereas  the  neutrino  energy,  the  baseline  length,  and  the  matter- 
induced  effective  potential  vary  from  experiment  to  experiment.  The  present  best-£t  values 
and  3(7  allowed  ranges  of  the  fundamental  neutrino  parameters  found  in  a  recent  global  fit 
of  the  neutrino  oscillation  data  [40]  are  summarized  in  Table  1.  Unless  otherwise  stated, 
all  calculations  in  the  present  paper  are  performed  for  the  following  values  of  the  neutrino 
parameters:  =  7  ■  10“^  eV^,  ^12  =  33°,  and  623  =  45°.  For  the  atmospheric  mass 

squared  difference  Amg^  we  adopt  the  current  best-fit  value  given  by  the  Super-Kamiokande 
Collaboration,  |Am|g|  =  2  ■  10“^  eV^  [2],  which  is  slightly  smaller  than  the  value  given  in  Ta¬ 
ble  1.  The  sign  of  Am|g  is  related  to  the  neutrino  mass  hierarchy:  for  the  normal  (inverted) 
hierarchy  one  has  Am|]^  >  0  (A?7i|g  <  0).  For  the  leptonic  mixing  angle  ^13  we  allow  values 
below  the  90%  C.L.  upper  bound  found  in  the  global  fit  of  the  neutrino  oscillation  data  [40] 
for  lAmgj^l  fixed  at  2  ■  10“^  eV^: 

6'i3<10.8°,  or  si3  =  sin  6^13  <  0.19,  or  sin^  26^13  <  0.14 .  (3) 

For  the  CP  violation  phase  (5cp  we  allow  values  between  0  and  27r. 
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Parameter 

Best-£t  value 

Range  (3a) 

Am^^ 

6.9-  10-5  eV^ 

(5.4  d- 9.5)  ■  10-5  eV^ 

[Am^il 

2.6  ■  10-3  eV^ 

(1.5  d- 3.7)  ■  10-3  eV^ 

^12 

33.2° 

28.6°  d- 38.6° 

013 

4.4° 

Od- 13.4° 

to 

CO 

46.1° 

33.8°  d- 58.1° 

5cp 

- 

0d-27r 

Table  1:  Present  best-fit  values  and  3a  allowed  ranges  of  the  fundamental  neutrino  param¬ 
eters  from  a  three- flavor  fit  to  global  neutrino  oscillation  data  [40]. 


Inspecting  the  values  of  the  fundamental  neutrino  parameters  in  Table  1,  one  can  identify 
two  natural  candidates  for  small  expansion  parameters  of  the  neutrino  oscillation  probabil¬ 
ities.  These  are  the  small  leptonic  mixing  angle  613  (or,  equivalently,  S13)  and  the  mass 
hierarchy  parameter 


Amii 


~  0.026 , 


(0.018)  0.021  <  a  <  0.036  (0.053)  at  90%  C.L.  (3a) , 


(4) 


where  we  have  taken  the  current  best-fit  value  and  allowed  ranges  from  Ref.  [40].  In  the 
following,  we  will  derive  a  number  of  formulas  for  series  expansions  of  the  neutrino  oscillation 
probabilities  in  these  small  quantities.  Comparing  Eqs.  (3)  and  (4),  one  realizes  that  current 
data  constrain  the  parameter  a  to  a  relatively  narrow  range,  while  S13  is  only  bounded  from 
above  and  might  be  signihcantly  larger  or  smaller  than  a.  The  relative  size  of  these  two 
expansion  parameters  will  be  important  for  the  validity  of  a  given  type  of  expansion. 

In  order  to  hnd  the  neutrino  oscillation  probabilities  for  a  given  experimental  setup,  one 
has,  in  general,  to  solve  the  Schrodinger  equation  for  the  neutrino  vector  of  state  in  the  flavor 
basis  \u{t))  =  (ueit)  u^(t)  z/^(t))'^: 


with  the  effective  Hamiltonian 


z/(t))  =  H\u{t)) 


(5) 


H  ~ 


1 


U  diag(0,  Am2i,  +  diag(l/',  0,  0) . 


(6) 


Here  V  is  the  charged-current  contribution  to  the  matter-induced  effective  potential  of  i/g  [25]. 
We  have  disregarded  the  neutral-current  contributions  to  the  neutrino  potentials  in  matter, 
since  they  are  the  same  for  i/g,  z/^,  and  i/r  ^  and  so  do  not  affect  neutrino  oscillations.  Note  that 
Eq.  (6)  holds  for  neutrinos,  whereas  for  antineutrinos  one  has  to  perform  the  replacements 


u  ,  V^-V. 


(7) 


^Up  to  tiny  radiative  corrections  [41]  which  are  negligible  except  at  very  high  densities  available,  e.g., 
inside  supernovae. 
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The  potential  V  (x)  is  given  in  convenient  units  by 

V(x)  ~  7.56  X  10“^^  f  Ve(x)  eV  ,  (8) 

\g/cm-^  / 

where  p(x)  is  the  matter  density  along  the  neutrino  path  and  Velx)  is  the  number  of  electrons 
per  nucleon.  For  the  matter  of  the  Earth  one  has,  to  a  very  good  accuracy,  17  —  0.5. 

For  many  practical  applications  (such  as  long-baseline  accelerator  experiments,  as  well 
as  oscillations  of  atmospheric,  solar,  and  supernova  neutrinos  inside  the  Earth  when  they 
do  not  cross  the  Earth’s  core)  it  is  a  very  good  approximation  to  assume  that  the  matter 
density  along  the  neutrino  trajectory  is  constant  (see,  e.g.,  Refs.  [42-44]).  Typical  values  for 
the  matter  density  are  pcmst  —  3g/cm^  in  the  Earth’s  crust  and  Pmantie  —  4.5g/cm^  in  its 
mantle.  In  situations  where  the  neutrinos  also  cross  the  Earth’s  core  or  for  strongly  varying 
matter  density  profiles  like  those  inside  the  Sun  or  supernovae,  the  constant  matter  density 
approximation  is  not  valid. 

The  neutrino  oscillation  probabilities  can  be  found  as  Fag  =  |S'pa(f,  to)  P,  where  S(t,to) 
is  the  evolution  matrix  such  that 

|z/(t))  =  S'(t,to)k(to)) ,  S{to,to)  =  l.  (9) 

Note  that  S(t,to)  satisfies  the  same  Schrodinger  equation,  Eq.  (5),  as  \i^(t)).  In  the  case 
of  matter  of  constant  density,  the  evolution  matrix  can  be  obtained  by  diagonalizing  the 
Hamiltonian  in  Eq.  (6)  according  to  H  =  U'HU'\  where  U'  is  the  leptonic  mixing  matrix  in 
matter  and  H  =  diag(i7i,  7^2,  E^).  The  evolution  matrix  is  then  given  by 

3 

Spait,to)  =  ,  a,/3  =  e,/i,r,  (10) 

i=l 

where  we  have  identified  L  =  t  — 

Before  presenting  our  results  for  the  neutrino  oscillation  probabilities,  we  discuss  some  of 
their  general  properties  as  well  as  relations  between  them.  First,  we  note  that  Eq.  (7)  implies 
that  one  can  relate  the  oscillation  probabilities  for  antineutrinos  to  those  for  neutrinos  by 

Pag  =  Pafsi^cp  ^ -Scp,  V  ^ -V) ,  a,(3  =  e,p,T.  (11) 

Second,  in  general,  i.e.,  both  in  vacuum  and  in  matter  with  an  arbitrary  density  profile,  it 
follows  from  the  unitarity  of  S'(f,to)  (conservation  of  probability)  that 

^  ^  Pag  ^  ^  Pag  1  )  Q.^  (3  6,  .  (12) 

a  g 

These  relations  imply  that  five  out  of  the  nine  neutrino  oscillation  probabilities  can  be  ex¬ 
pressed  in  terms  of  the  other  four  [45]. 

Besides  these  general  properties,  there  exists  an  additional  symmetry  due  to  the  specific 
parameterization  of  the  mixing  matrix  given  in  Eq.  (1)  and  the  fact  that  the  rotation  matrix 


5 


023  commutes  with  the  matter  potential  term  of  the  Hamiltonian  in  Eq.  (6).  It  is  easy  to 
show  that  the  evolution  matrix  can  be  written  as 


S {t,  to)  =  023S\t,  to)023  (13) 

where  S'(t,to)  does  not  depend  on  623.  This  can  be  used  to  prove  some  useful  relations 
between  the  probabilities.  Let  us  denote  ^ 

^a/3  =  Th/3 (^23  ^  C23  ,  sin 26^23  ^ -sin 26^23 ) ,  a,(3  =  e,n,T.  (14) 

Using  Eqs.  (13)  and  (14),  one  can  readily  show  that 

Per  =  Pefi  ,  Pth  =  Pfir  ,  Ptt  =  Pfifi  ,  (15) 


while  Pee  turns  out  to  be  independent  of  6*23. 

Out  of  the  three  conditions  in  Eq.  (15),  only  two  are  independent,  as  each  of  them 
can  be  derived  from  the  other  two  and  the  unitarity  conditions  (12).  Hence,  the  number 
of  independent  neutrino  oscillation  probabilities  is  reduced  to  two.  Thus,  we  come  to  the 
important  conclusion  that  all  the  nine  neutrino  oscillation  probabilities  can  be  expressed 
through  just  two  independent  probabilities  provided  that  their  dependence  on  the  mixing 
angle  623  is  known.  However,  the  choice  of  these  independent  probabilities  is  restricted: 
they  should  not  include  Pee,  which  is  independent  of  6*23;  nor  should  they  be  a  pair  of  the 
probabilities,  which  are  T-reverse  of  each  other,  or  go  into  each  other  (or  T-reverse  of  each 
other)  nnder  the  transformation  S23  0^3  ,  sin  2023  — ^  —  sin  2023- 

One  possible  choice,  and  the  one  that  we  will  use,  is  Pe^  and  P^r-  For  completeness,  we 
give  here  the  expressions  for  all  the  other  neutrino  oscillation  probabilities  in  terms  of  these 
two.  Using  Eqs.  (12)  and  (15),  one  easily  hnds 


p 

^  ee 

—  1  —  {Peii  +  PeiJt)  , 

(16) 

P 

^  er 

=  Pe^L  , 

(17) 

P 

^  /te 

h 

+ 

1 

(18) 

P 

MM 

=  1  -  P  -  P 

(19) 

P 

^  re 

=  P  +  P  -  P 

(20) 

P 

^  TfJ, 

=  P^MT  , 

(21) 

P 

1  7-7- 

=  1  —  {Pen  +  Pfir)  ■ 

(22) 

In  addition  to  the  above  relations,  one  can  stndy  the  transformations  of  the  neutrino 
oscillation  probabilities  nnder  the  time  reversal  Pa/s  — >  Ppa-  It  can  be  shown  [28]  that  in 
matter  with  an  arbitrary  density  prohle 

=  Th/3(5cp  ^ -(5cp,  U(a;)  ^  Urev(a;)) ,  a,P  =  e,fi,T,  (23) 

^The  transformation  in  Eq.  (14)  can  be  achieved,  e.g.,  through  the  substitution  623  O23  +  7r/2  or 

023  ^  023  +  37r/2. 
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where  Kev(a^)  is  the  “reverse”  potential,  which  corresponds  to  the  interchanged  positions 
of  the  neutrino  source  and  the  detector.  In  the  case  of  symmetric  matter  density  prohles 
(including  matter  of  constant  density),  V^ew{.x)  =  V{x),  and  Eq.  (23)  simplihes  to 

Ppa  =  Pa(}{.^cv  ^ ,  a,/3  =  e,fi,T.  (24) 

While  Eq.  (24)  does  not  further  reduce  the  number  of  independent  probabilities,  it  yields 
relations  which  can  be  useful  for  cross-checking  the  formulas  for  in  the  case  of  matter 
with  symmetric  density  prohles. 

By  applying  the  rule  given  in  Eq.  (11),  one  can  obtain  from  Eqs.  (16)-(22)  the  cor¬ 
responding  probabilities  for  the  antineutrino  oscillations.  Thus,  the  expressions  for  all  18 
probabilities  of  neutrino  and  antineutrino  oscillations  can  be  found  from  the  formulas  for 
just  two  independent  neutrino  oscillation  probabilities,  which,  as  was  already  mentioned,  we 
choose  to  be  Pe^  and  P^r-  In  order  to  be  more  explicit,  we  will  in  some  cases  also  give 
formulas  for  additional  neutrino  oscillation  channels. 

In  the  following  Secs.  3,  4,  and  5,  we  give  our  results  for  the  various  series  expansions  of 
the  neutrino  oscillation  probabilities.  We  will  adopt  the  following  abbreviations: 

(25) 

(26) 

3  Series  expansion  up  to  second  order  in  o:  and  S13 


^  Amoi  L 
A  = 


A  = 


4E  ’ 
2EV  VL 


Amii 


3.1  Matter  of  constant  density 


In  this  section,  we  present  the  series  expansion  formulas  for  three-havor  neutrino  oscilla¬ 
tion  probabilities  in  matter  of  constant  density  up  to  second  order  in  both  a  and  S13.  The 
probabilities  are  calculated  by  diagonalizing  the  Hamiltonian  (6)  up  to  second  order  in  these 
parameters,  as  described  in  Appendix  A.l.  We  hnd  for  the  eigenvalues  of  the  Hamiltonian 


El  ~ 
E2  ~ 

P3  ~ 


^^31 

2E 

2E 

2E 


A  -|-  Q;  s-in  -|-  s 


A  -^-2 


2 


-|-  O 


sin"^  2012 
4A 


Of  0^2  ^ 


sin  26*12 


4A 


1  —  s 


A 


13 


A  -  1 


(27) 

(28) 
(29) 


Calculating  the  mixing  matrix  in  matter  U'  as  described  in  Appendix  A.l  and  using  Eq.  (10) 
for  the  evolution  matrix  S',  it  is  straightforward  (although  somewhat  tedious)  to  derive  the 
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following  expressions  for  the  neutrino  oscillation  probabilities: 


Ppp  =  1  —  a  sin  29 


sinMA  2 


-4s 


(A -1)2  ’ 


n  2  •  2ri/i  2  sin^AA  2  2  siii^fA  -  1)A 
Pef,  =  a  sin  2^12  C23  +  4  s^g  ^23  _  ^^2 


2  a  Si3  sin  26*12  sin  26*23  cos(A  —  (5cp) 


sinAA  sin(A  —  1)A 
“4  A-  1 


sin^  AA  .22  sin^(A  —  1)A 


Per  =  a  sin  2012  S23 - -  +  4  Si3  C23 


2  a  S13  sin  20i2  sin  2023  cos(A  —  (5cp) 


(A -1)2 

sin  AA  sin(A  —  1)A 
“4  A^l 


Puu  =  1  —  sin^  2023  sin^  A  +  a  c?,  sin^  2023  A  sin  2  A 


sin^  AA 


—  a  sin  2012  C23 - -7^ - a  C12  sin  2023  A  cos  2  A 


1  9  9  ^  9  ^  /  *  sin  A  A  ,  ,  ,  .  A 

H - -  sin^  20i2  sin^  2023  (  sin  A  — - —  cos(A  —  1)  A - sin  2A 

2A  \  A  2 


_  2  2  sin  (A-l)A 

4^13^23 

- —  s?3  sin^  2023  ( sin  A  cos  AA  — 11^  _  sin  2 A 

A  -  1  M  A  -  1  2 


2  a  S13  sin  20i2  sin  2023  cos  (5cp  cos  A 


sin  AA  sin(A  —  1)A 
“4  A-  1 


2  ^  sin  /\  'X 

+  — — ^q;Si3  sin  2012  sin  2023  cos  2023  cos^cp  sin  A  (A  sin  A - — —  cos(A  —  1)A  j  , 


Par  =  sin^  2023  sin^  A  —  a  c?o  sin^  2023  A  sin  2  A  +  cfn  sin^  2023  A^  cos  2  A 


- -  o?  sin^  2012  sin^  2023  (  sin  A 


sin  A  A  ,  ^  ^ .  A  .  ^  . 

— - —  cos(A  —  1)A - sin2A 

A  2 


2  ,  j  a  (  A  A  A  sin(A  —  1)A  A  .  . 

+  — - s?o  sin^  2023  I  sin  A  cos  AA - - - A  sin  2 A 

A-1^^  M  A-1  2 


2q;si3  sin  2012  sin  2023  sin(5cp  sin  A 


sin  AA  sin(A  —  1)A 
“4  A-1 


— — sin  26*12  sin  26*23  cos  26*23  cos(5cp  sin  A  (  A  sin  A 


sin  A  A 


(34) 


Formally,  our  calculations  are  based  upon  the  approximations  a,  S13  1  and  no  explicit 

assumptions  about  the  values  oi  L/ E  are  made.  However,  we  remark  that  the  series  expansion 
formulas  (30)-(34)  are  no  longer  valid  as  soon  as  aA  =  Arri^iL/^AE)  becomes  of  order  unity, 
i.e.,  when  the  oscillatory  behavior  due  to  the  “solar”  mass  squared  difference  Am^i  becomes 
relevant.  This  can  happen  for  very  long  baselines  and/or  very  low  energies.  See  also  Sec.  6 
for  a  detailed  discussion  of  the  accuracy  of  these  formulas. 

From  Eqs.  (27)-(29)  one  can  see  that  in  vacuum  {A  =  0)  and  at  the  atmospheric  reso¬ 
nance  [A  =  1)  the  expressions  for  the  eigenvalues  are  divergent,  and  one  would  expect  the 
expansion  to  break  down.  In  these  cases,  two  out  of  the  three  eigenvalues  of  the  unperturbed 
Hamiltonian  are  degenerate  and,  strictly  speaking,  the  degenerate  perturbation  theory  rather 
than  the  ordinary  one  should  be  employed.  However,  it  turns  out  that,  though  the  eigenvalues 
(27)-(29)  are  divergent,  the  neutrino  oscillation  probabilities  are  hnite  in  the  limits  A  — 0 
and  A  — 1.  The  reason  for  this  interesting  behavior  is  a  cancellation  of  divergences  between 
the  eigenvalues  and  the  matrix  elements  of  the  leptonic  mixing  matrix  in  the  calculation  of 
the  evolution  matrix  according  to  Eq.  (10).  In  particular,  in  the  limit  A  — 0,  the  correct 
vacuum  neutrino  oscillation  probabilities  are  obtained. 

We  shall  now  compare  the  above  results  with  those  existing  in  the  literature.  Equa¬ 
tions  (27)-(32)  have  previously  been  derived  in  Ref.  [35]  and  conhrmed  in  Ref.  [36].  Expres¬ 
sions  (33)  and  (34)  are  new.  In  Ref.  [38],  an  expression  for  was  found  to  hrst  order  in 
a,  including  the  0{asi3)  term,  which  can  be  compared  with  the  corresponding  terms  in  our 
Eq.  (34).  We  hnd  that,  while  our  0{1)  and  0{a)  terms  coincide  with  those  in  Eq.  (A3)  of 
Ref.  [38],  our  0{asi^)  term  is  quite  different.  In  particular,  we  disagree  with  the  statement 
in  Ref.  [38]  that  to  order  a  the  probability  P^r  does  not  depend  on  the  CP-violating  phase 
(5cp.  The  existence  of  a  term  proportional  to  asin^cp  in  P^r  is  actually  expected,  since  in 
matter  of  constant  density  the  phase  (5cp  is  the  sole  source  of  T-violation,  and  up  to  the  sign, 
the  T-odd  terms  in  all  three  transition  probabilities  must  be  the  same  [26].  This  is  indeed 
seen  in  Eqs.  (31),  (32),  and  (34).  We  also  note  that  the  term  of  order  asis  in  Eq.  (A3)  of 
Ref.  [38]  diverges  at  the  atmospheric  resonance  (A  =  1),  while  our  expression  (34)  is  regular 
at  all  physical  values  of  parameters. 

3.2  Vacuum  neutrino  oscillation  probabilities  up  to  second  order  in  a  and  S13 

The  vacuum  neutrino  oscillation  probabilities  up  to  second  order  in  a  and  S13  are  given  by 
Pl^  =  1-0^  sin^  2^12  A^  -  4  sin^  A  ,  (35) 

sin^  26*12  C23  A^  -|-  4  S13  S23  sin^  A  -|-  2  a  S13  sin  26*12  sin  26*23  cos(A  —  (5cp)  A  sin  A  , 

(36) 

Pg™'^  =  sin^  26*12  -S23  A^  -|-  4  C23  sin^  A  —  2  a  S13  sin  26*12  sin  26*23  cos  (A  —  (5cp)  A  sin  A  , 

(37) 
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=  1  -  sin^  2023  sin^  A  +  acj^  sin^  2023  A  sin 2 A 

—  [sin^  20i2  cl^  +  c^2  2023  (cos  2 A  —  5^2)]  +  4  s^g  S23  cos  2023  sin^  A 

—  2  a  Si3  sin  20i2  S23  sin  2023  cos  (5cp  A  sin  2  A  ,  (38) 

=  sin^  2023  sin^  A  —  a  0^2  sin^  2023  A  sin  2  A 
+  sin^  2023  A^  ^0^2  cos  2  A  —  ^  sin^  20i2  sin^  A^ 

—  2  Sgg  sin^  2023  sin^  A 

+  2  a  Si3  sin  20i2  sin  2023  (sin  (5cp  sin  A  —  cos  2023  cos  (5cp  cos  A)  A  sin  A  .  (39) 

4  Series  expansion  up  to  first  order  in  ol 

In  this  section,  we  expand  the  probabilities  up  to  first  order  in  the  small  parameter  a  while 
keeping  their  exact  dependence  on  0i3.  One  expects  these  formulas  to  be  useful  for  relatively 
large  values  of  S13.  In  addition,  they  will  correctly  account  for  the  “atmospheric”  resonance 
driven  by  the  parameters  Am|g  and  0i3. 

4.1  Matter  of  constant  density 

The  eigenvalues  of  the  Hamiltonian  (6)  up  to  first  order  in  a  are 

(40) 

/\yy)  1 

-^2  —  2  ~  ^13)  +  ^^^Q('Si2  (Cl3  +  1  —  AcOs20i3)  ,  (41) 

/\fY}  1 

E3  —  -  (1  +  A  +  O13)  +  a  s‘^2  (Ci3  —  1  +  A  cos  2013)  ,  (42) 

where 

C*i3  =  \J sin^  2013  +  (A  —  cos  2013)^  .  (43) 

To  calculate  the  probabilities  to  first  order  in  a  we  followed  two  different  approaches:  First, 
we  used  the  Cayley-Hamilton  formalism  as  described  in  Appendix  A. 2,  to  series  expand  the 
evolution  matrix,  and  second,  we  considered  the  const  ant- density  limit  of  the  expansion  of 
the  evolution  equation  described  in  Appendix  B.l.  Both  methods  gave  the  same  results, 
which  is  a  useful  cross-check  of  our  calculations.  Writing 

=  +  +  (44) 
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we 


obtain  the  following  expressions  for  the  channel; 


^13 


p(i)  o  2  sin2  20i3  . 

PeJ  =  2si2  —pT2 - Sin  ^13  A 

O' 13 


X 


A 


13 

cos  Ci3  A 


Similarly,  for  the  v^ 


/I  0/3  ^  sinC'i3A  COS2013  -  A' 
(1  -  A  cos  2013)  -  A — - - - - 

^13  ^13 


e  ^ 


Ci3 

channel  we  hnd 


(0)  2  sin2  20i3  .  2n  \ 

=  ^23  —^2 -  sin  C'i3A  , 

0'13 

CO  9  9  sin^20i3  ^  . 

PlJ  =  -2512^23  —7^2 -  SinC'i3A 

0'13 


X 


A 


cosCisA  sin  Ci3  A  cos  2013  ■ 

- - (1  -  A cos 2013)  -  A— - ^ 

0'13  0/13  0-13 


2013  -  A' 


I  sin  (5cp  [cos  C13A  -  cos(l  +  A)  A]  C13 


sin 

+  Si3  sin  2012  sin  2023 — ^^2 
+  cos  (5cp  [Ci3  sin(l  +  A)  A  —  (1  —  A  cos  20i3)  sin  C13  A]  | , 


and  hnally,  for  the  Uj.  channel  we  have 


pW  =  -sin2  2023 


1  - 


cos  2013  —  A 


C 


13 


-f  siri'^  1(1  +  /I  -  Ci3)A 


/  3082/^13-/11  .  ol/,  ,  ^  lA  lsill^2fl 

(1  +  j  +  -^  +  -  2  ^ 

dt  ~  ~2  a|2  C12  -  812813^^(1  +  28^3 A  +  A‘^) 


lsm^26li3  .r  ' 

-^85 — sm  C13A 


(45) 


(46) 


(47) 


(48) 


(49) 


COS  C13A  sin(l  +  A)  A 


+  2  [c?2Ci3  -  cl^sls  +  s?2«?3  +  “  c^)  cos(l  +  A)A 

0-13 

2  sin^  2013  sin  C13  A 
+  ^12  — ^^2 - 


C 


13 


X 


13 

A  .  ,A  ^  Acos20i3  — A  .  ^  .  A,  3  3,31  A  ^  A 

—  sin(l  +  A)  A  +  — - — - sin  PisA  —  (1  —  A  cos  20i3)  cos  PisA 

A  A  0/13 

—  [cosCisA  —  cos(l  +  A)  A] 

Ai3 

2 


Si3  sin  2012  sin  2023  f  9  ,,  sinC 

+  2/li  "|2c;38.n4eP— 


Z21Ci3  f  0-13 

—  cos 2023  cos(5cp(l  +  A)  [cosCisA  —  cos(l  +  A)A]^ 
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+  cos  26*23  cos  ^CP 


•  ^  \  A  cos  26*13  —  ^  ^  A 

sin(l  +  A) A  H - - sin(7i3A 


C 


13 


X 


sin  Cl  3  A 


(1  +  24a  +  -  (1  +  A)  sin(l  +  A) A 


C 


13 


(50) 


In  the  above  formulas,  one  may  identify  the  effective  mixing  angle  in  matter  in  the  1-3  sector 
6*'^3,  which  is  determined  through  the  expressions 


sin  26*^3  = 


sin  26*13 


cos26*'^3  = 


cos  26*13  ~  ^ 


C 


13 


(51) 


appearing  frequently  in  Eqs.  (45)-(50).  Furthermore,  the  combination  C13A  appearing  as  the 
argument  of  sine  or  cosine  corresponds  to  the  effective  A  in  matter.  In  the  limit  when  ^13  is 
small,  one  has  C13  ~  A  —  1,  and  expanding  Eqs.  (45)-(50)  up  to  second  order  in  S13  yields 
the  double  expansions  given  in  Sec.  3,  except  for  the  terms  of  order 

Equations  (40)-(42)  and  (45)-(48)  have  previously  been  derived  in  Ref.  [35]  and  conhrmed 
in  Ref.  [36].  Expressions  (49)  and  (50)  are  new. 


4.2  Vacuum  oscillation  probabilities  up  to  first  order  in  ex 


Taking  the  limit  A  — 0  in  Eqs.  (45)-(50),  it  is  straightforward  to  obtain  the  neutrino  oscil¬ 


lation  probabilities  in  vacuum  to  first  order  in  a: 

pjo)vac  ^  ^  20^^  sin2  A  ,  (52) 

pa)vac  _  ^^2^  gi^2  20^^  sin  2  A  ,  (53) 

P4™^=sin2  20i34sin2A,  (54) 

=  -A4  sin^  20134  sin  2A 

-|-  A  sin  2012S134  sin  2023  (2  sin  5cp  sin^  A  -|-  cos  (5cp  sin  2  A) ,  (55) 

pWvac  _  ^4^  gi^2  20^^  sin2  A  ,  (56) 

pa)vac  ^  -Acl^  sin^  2023(4  -  sin2A 

-|-  A  sin  2012  'Si34  ^023  (2  sin  (5cp  sin^  A  —  cos  (5cp  cos  2023  sin  2  A) .  (57) 


5  Series  expansion  up  to  first  order  in  S13 

In  this  section,  we  expand  the  probabilities  up  to  first  order  in  the  small  parameter  S13  while 
keeping  their  exact  dependence  on  a.  These  formulas  are  expected  to  be  useful  whenever 
neutrino  oscillations  driven  by  the  solar  mass  squared  difference  Am^i  are  important.  In 
practice,  this  means  low  neutrino  energies  and  long  baselines. 
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5.1  Matter  of  constant  density 


The  eigenvalues  of  the  Hamiltonian  (6)  up  to  first  order  in  S13  are  given  by 


^  Amoi 
El  ~ - ^ 

2E 

Ami 


Eo  ~ 


31 


2E 


^  +  2(1 -^12) 

^  Cl/,  ^  X 

-  +  -(1  +  C^12) 


E;,  ~ 


^™ll 

2E 


(58) 

(59) 

(60) 


where 


C 


12 


^cos  29 


2 


(61) 


Note  that  these  eigenvalues  are  independent  of  the  expansion  parameter  S13,  which  is  consis¬ 
tent  with  Eqs.  (27)-(29),  where  the  lowest-order  in  S13  corrections  to  the  eigenvalues  appear 
only  at  order  5^3. 

As  in  the  case  of  the  single  expansion  in  a,  we  have  calculated  the  probabilities  using 
two  different  methods:  by  series  expanding  the  evolution  matrix  using  the  Cayley-Hamilton 
formalism  as  described  in  Appendix  A. 2,  and  using  the  constant-density  limit  of  the  single 
expansion  in  S13  of  the  evolution  equation  described  in  Appendix  B.2.  We  have  found  identical 
results.  Writing  the  probabilities  as 


=  d?  +  S.3  di'  +  ,  (62) 

we  obtain  the  following  expressions  for  the  z/g  — i^e  channel; 

Pee^  =  1  -  aCuA  ,  (63) 

^12 

Pi‘'=0.  (64) 


The  absence  of  any  Erst  order  corrections  to  Pge  is  consistent  with  Eq.  (30).  Similarly,  for 
the  z/g  — z/^  channel  we  find 


(0)  ^  2  sin^  2^12 

■^e/4  *^23  ^2 

'-'12 


sin^  q;Ci2  A , 


(65) 


(3)  sin  26*12  .  (1  -  a)  sinaCi2A 

P},/  =  — - sm  26*23 


efi 


c 


12 


1  —  A  —  a  -|-  Aacy 


12 


—  cos  (5cp 


sin(A  +  a  —  2)  A  —  sin  q;Ci2  A 


sin  (5cp  [cos  q;Ci2 A  —  cos(A  +  a  —  2) A] 
cos  26*12  -  ^  q;ACi2  sin^  26* 


712 


C 


12 


2(1  -  a)  Cf, 


(66) 
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and  for  the 

p{0) 

At-r 


p{l) 

^lT 


Vt-  channel  we  have 


1 

2 


=  -  sin^  29 


23 


1  -  - 


1  sin^  29 


12  .  2 


2  C^2 


sin^  aCi2^  —  cos(q;Ci2  +  A  +  a  —  2)  A 


,  cos  26*12  —  -  ,  ,  / 

—  I  1 - ^  1  sin  q;Ci2A  sin(A  +  a  —  2)  A 


C 


12 


sin  2^12  . 

sin  26*23 


X 


1  —  A  —  a  +  Aacf2 

r 

cos  26*23  cos  (5cp 


X 


<^12 

I"  q;AC*i2 

cos  26*12  -  f 
<^12 

cos  26*12  -  i  ,  2(1  -  a) 


C 


(cosaCuA  —  cos(A  +  a  —  2)  Ay 
sinQ;Ci2A  +  sin(A  +  a  —  2)  A 

/ 

sin  q;Ci2A  +  sin(A  +  a  —  2)  A 


+ 


12 


aAC 


12 


+  sin  (5cp  (1  —  a)  (cos  q;Ci2A  —  cos(A  +  a  —  2)  A)  sin  q;Ci2A  >  . 


(67) 


(68) 


In  this  case,  one  may  identify  in  Eqs.  (63)-(68)  the  effective  “solar”  mixing  angle  in  matter 
6*^2,  which  is  determined  by 


sin  26*^2  = 


sin  20 


12 


c 


COs20'^2  = 


cos  2012  — 


12 


c 


12 


(69) 


The  combination  q;C'i2A  appearing  as  argument  of  sine  or  cosine  corresponds  to  oscillations 
with  the  “solar”  frequency  in  matter.  Furthermore,  we  note  that  in  the  limit  when  a  is  small, 
one  has  C12  A/ a  —  cos20i2,  and  expanding  Eqs.  (63)-(68)  up  to  second  order  in  a  yields 
the  double  expansions  given  in  Sec.  3,  except  for  the  terms  of  order 

Neutrino  oscillations  probabilities  in  matter  of  constant  density  expanded  to  first  order 
in  Si3,  but  exact  in  a,  presented  in  this  subsection,  have  not  been  previously  published  and 
are  entirely  new. 


5.2  Vacuum  oscillation  probabilities  up  to  first  order  in  S13 

Taking  the  limit  A  — 0  in  Eqs.  (63)-(68),  one  recovers  the  vacuum  probabilities  expanded 
to  hrst  order  in  S13: 

P™’”  =  1  -  sin=  2»i2  sin^  aA , 

Piy“  =  0, 

pmvac  ^  si„2  _ 

—  cosiAcp  sill  2(?i2  sin  2(?23  [sin^  A  —  sin^(l  —  a)A  +  cos  20,2  sin^  aA] 


14 


(70) 

(71) 

(72) 


(73) 

(74) 


+  -  sin  (5cp  sin  26*i2  sin  2023  [—  sin  2 A  +  sin  2(1  —  a)  A  +  sin  2q; A] , 

26*23  sin^  A  +  c^2  sin^  26*23 [sin^(l  ~  ~  -5^2  , 

=  cos  (5cp  sin  26*12  sin  26*23  cos  26*23  [— sin^  A  +  sin^(l  —  q;)A  +  cos  26*12  sin^  a  A] 

+  -  sin  Scp  sin  26*12  sin  26*23  [~  sin  2 A  +  sin  2(1  —  a)  A  +  sin  2q; A] .  (75) 

6  Qualitative  discussion  and  tests  of  accuracy 

In  this  section,  we  discuss  the  qualitative  behavior  of  the  neutrino  oscillation  probabilities 
and  the  relevance  of  matter  effects  (Sec.  6.1).  We  also  assess  quantitatively  the  accuracy 
of  the  various  expansions  and  compare  each  type  of  expansion  with  an  exact  numerical  cal¬ 
culation  of  the  corresponding  probability  (Secs.  6.2  and  6.3).  We  examine  in  detail  which 
type  of  expansion  is  most  accurate,  depending  on  the  values  of  the  fundamental  parameters 
(especially  a  and  S13),  and  on  the  experimental  conhguration  characterized  by  the  neutrino 
energy  E  and  the  baseline  length  L.  In  Sec.  6.4,  we  discuss  issues  related  to  the  applica¬ 
tion  of  our  formulas  to  neutrino  oscillation  experiments.  As  we  have  shown  in  Sec.  2,  the 
probabilities  for  all  oscillation  channels  can  be  obtained  from  the  expressions  for  and  P^r 
[see  Eqs.  (16)-(22)].  Therefore,  we  focus  in  the  following  on  the  Ve  channel;  we  briefly 

comment  also  on  the  other  neutrino  oscillation  channels. 

6.1  The  relevance  of  matter  effects  and  the  probability 

Before  considering  the  accuracy  of  the  expansion  formulas,  we  discuss  in  this  subsection 
the  features  and  the  relevance  of  matter  effects  for  the  oscillation  probability  P^fj,. 

The  discussion  applies  also  to  Per]  matter  effects  on  the  probabilities  Prr,  and  P^r  are 
signihcantly  weaker  than  those  on  P^^  and  Per-  Figure  1  shows  contours  of  Pe^  for  matter 
of  constant  density  calculated  numerically  from  Eq.  (6)  without  any  approximations,  for  a 
wide  range  of  baseline  lengths  and  neutrino  energies.  Many  features  of  this  figure  can  be 
understood  by  considering  the  expression  for  the  two-flavor  neutrino  oscillation  probability 
in  matter  of  constant  density: 

sin^  20  t 

P  =  sin^  C^2u ,  P  =  ysin^  26* -f  (cos26*  -  A2i.)2  ,  (76) 

where,  in  analogy  to  Eqs.  (25)  and  (26),  we  dehne  /S.2y  =  Arn^L/{4E),  A2U  =  2EV/Am‘^, 
and  6*  and  Am^  are  the  generic  two-flavor  neutrino  oscillation  parameters.  From  Eq.  (76) 
it  is  clear  that  for  A2U  ^  1  one  has  P  ~  1,  and  vacuum  neutrino  oscillations  with  P  ~ 
sin^  26*  sin^  A21,  are  recovered.  For  A2U  =  cos  26*  one  can  see  in  Eq.  (76)  also  the  MSW 
resonance  [25],  which  leads  to  P  =  sin 20  in  Eq.  (76)  and  to  the  effective  mixing  angle  in 
matter  sin^  20'  =  sin^  20/P^  =  1.  The  general  resonance  conditions  for  three  flavors  are  much 
more  complicated  than  those  in  the  two-flavor  case.  However,  due  to  the  hierarchy  of  the 
mass  squared  differences  and  smallness  of  0i3,  one  can  use,  to  a  very  good  approximation,  the 
respective  effective  two-flavor  picture.  This  leads  to  the  “solar”  and  “atmospheric”  resonance 
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Figure  1:  Contours  of  Pg/i  calculated  numerically  without  approximations  for  sin^  2^13  = 
0.02,  =  2  ■  10“^  eV^,  a  =  0.03,  6*12  =  33°,  6*23  =  45°,  and  5cp  =  0.  Matter  of  constant 

density  p  =  3g/cm^  is  assumed,  and  the  probability  is  averaged  over  a  Gaussian  energy 
resolution  of  1%. 


energies  shown  as  vertical  lines  in  Fig.  1 


A 

a 


2EV 

/Aml^ 


cos  2012  , 


2EV 

Amii 


cos  2013  • 


(77) 


The  “atmospheric”  resonance  is  clearly  visible  in  Fig.  1  at  Pres  —  8.6  GeV  and  L  >  5500  km. 
The  “solar”  resonance  occurs  at  the  energy  Pres  —  0.11  GeV.  Note  that  the  solar  resonance  is 
not  as  pronounced  in  Fig.  1  as  the  atmospheric  one,  since  the  neutrino  oscillation  amplitude 
in  vacuum  sin^  20i2  is  already  quite  large. 

Far  enough  to  the  left  of  the  vertical  lines  in  Fig.  1  marking  the  solar  and  atmospheric 
resonance  energies  one  has  vacuum  neutrino  oscillations  with  the  “solar”  parameters  0i2, 
Am^i  and  with  the  “atmospheric”  parameters  0i3,  Arugi,  respectively.  Indeed,  the  typical 
LjE  pattern  is  clearly  visible  in  Fig.  1  to  the  left  of  the  vertical  lines.  The  diagonal  lines  in¬ 
dicate  the  constant  values  of  L/P  corresponding  to  the  hrst  solar  and  atmospheric  oscillation 
maxima  in  vacuum,  given  by  the  conditions 


aA 


Am^iT 

4P 


71 
2  ’ 


A 


4P 


71 
2  ’ 


(78) 


respectively.  To  the  right  of  the  vertical  lines  the  probability  is  dominated  by  matter  effects. 
For  A21,  3>  1  we  have  C  ~  A2U  in  Eq.  (76),  which  has  two  consequences:  First,  since 
A21/A21/  =  VL/2,  the  oscillation  frequency  becomes  independent  of  energy,  and  second,  the 
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oscillation  amplitude  becomes  suppressed  at  high  energies  because  sin^  29! oc  1/ E"^.  Both 
these  effects  are  apparent  in  Fig.  1  for  neutrino  oscillations  with  the  solar  as  well  as  the 
atmospheric  frequency. 

In  addition  to  the  dependence  on  the  neutrino  energy,  matter  effects  depend  crucially  on 
the  baseline.  From  Eq.  (76)  one  can  see  that  the  two-flavor  oscillation  probabilities  approach 
the  vacuum  ones  for  0^2^  n/2.  Far  above  the  MSW  resonance  energy,  i.e.  for  A2u  3>  1, 

one  hnds  C  ~  A2u-  Hence,  the  short-baseline  limit  is  CA2U  —  H21/A21/  =  HL/2  7r/2,  leading 
to  vanishing  matter  effects  for  L  L^aat  =  tt/H  ~  5453  km.  The  horizontal  lines  in  Fig.  1 
indicate  this  baselines  of  the  first  “matter  effect  maximum”  for  large  energies  at  L  =  Lmat., 
and  the  first  “matter  effect  minimum”,  L  =  2L^s.t  —  10907  km.  At  the  baseline  2Lmat;  where 
AA  =  HL/2  =  TT,  the  terms  proportional  to  A  and  asis  in  the  double  expansion  Eq.  (31) 
disappear,  and  only  the  term  proportional  to  survives.  Therefore,  this  baseline  is  especially 
useful  to  measure  S13,  since  ambiguities  due  to  parameter  degeneracies  are  avoided  [38].  A 
dedicated  analysis  of  this  “magic  baseline”  can  be  found  in  Ref.  [46].^  Note  that  the  distance 
-hmat  depends  only  on  the  matter  potential  V,  i.e.  the  energy  E  as  well  as  neutrino  masses 
and  mixings  do  not  enter. 

Let  us  now  consider  the  energy  region  well  below  the  MSW  resonance  energy,  where 
A2V  ^  1-  The  short-baseline  behavior  of  the  oscillation  probability  in  this  region  can  be 
understood  in  the  two-flavor  picture  by  expanding  Eq.  (76)  for  hxed  ^2^  {i.e.,  hxed  L/E) 
assuming  A2U  -C  1  and  A2uA2u  1-  One  hnds 

pmatter  _  pvacuum  ^  2  sin^  26  COS  26  siu  A2V  (siu  A2V  —  A2U  COS  A2u)  A2U  ■  (79) 

This  relation  shows  that  —  pvacuum  fie(;;reag0s  linearly  with  A2U  oc  E,  which  means 

that  the  matter  effects  vanish  for  very  small  energies,  as  mentioned  above.  In  addition,  for 
a  hxed  energy  the  right-hand  side  of  Eq.  (79)  becomes  zero  for  A2U  ^  1  or,  equivalently, 
L  Lvac  =  ^tiE/ Am?.  Therefore,  the  two  conditions  for  the  smallness  of  matter  ehects  due 
to  the  short  baseline  can  be  written  for  all  energies  as  L  Lmin  =  niin(Lvac,  Lmat)- 

The  relevance  of  matter  ehects  is  illustrated  in  Fig.  2.  In  the  upper  left  panel  we  show 
the  probability  diherence  for  the  case  a  =  0,  such  that  the  oscillations 

due  to  the  “solar”  frequency  are  switched  oh  and  only  matter  ehects  are  present  which  are 
related  to  the  “atmospheric”  resonance  energy  at  8.6  GeV.  It  can  be  seen  that  the  diherence, 
as  well  as  the  relative  diherence  (lower  left  panel),  vanish  for  L  Tmin.  The  left  plots  of 
Fig.  2  show  that  the  diherences  vanish  also  for  small  energies  and  long  baselines.  This  can 
be  understood  as  follows.  For  E  P^es  the  oscillation  amplitude  in  matter  is  always  close 
to  that  in  vacuum.  For  very  short  baselines,  this  is  also  true  for  the  oscillation  phases.  At 
intermediate  baselines,  the  matter-induced  correction  to  the  oscillation  phase,  though  much 
smaller  than  the  main  “vacuum”  term,  can  become  comparable  to  unity  and  so  cannot  be 
ignored.  However,  at  long  baselines  L  ^  Lvac  the  averaging  regime  sets  in,  which  implies  that 
any  shifts  in  the  oscillation  phase  become  unimportant  and  the  vacuum  oscillations  limit  is 
regained. 

^It  should  be  noted  that  the  quoted  number  for  L^at  holds  for  a  constant  matter  density  of  p  =  3g/cm^. 
It  will  differ  for  larger  matter  densities  like  those  in  the  mantle  of  the  Earth,  or  even  more  drastically 
in  astrophysical  or  cosmological  applications  of  neutrino  oscillations.  In  addition,  the  approximation  of 
constant  matter  density  is  then  frequently  not  justified.  For  a  realistic  Earth  matter  density  profile  one  finds 
2L„,at  ^  7250  km  [46]. 
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Figure  2:  Contours  of  -  Pe””l  panels),  and  -  PeT“''“l/^e"'' 

(lower  panels),  where  p™^aer  j-pvacuum^  jg  exact  oscillation  probability  in  matter  (vacuum) 

for  sin^26'i3  =  0.05,  =  2  ■  10“^  eV^,  and  623  =  45°.  In  the  left  panels  a  =  0,  whereas  in 

the  right  panels  a  =  0.03,  612  =  33°,  and  hcp  =  0.  Matter  of  constant  density  p  =  3g/cm^ 
is  assumed,  and  the  probabilities  are  averaged  over  a  Gaussian  energy  resolution  of  2%.  The 
straight  lines  indicate  the  relevant  resonance  energies,  hrst  oscillation  maxima,  and  the  hrst 
matter  effect  maximum  (see  also  Fig.  1). 


The  above  discussion  describes  the  regions  where  |  —  pvacuum  |  vanishes  or  becomes 

small  for  a  =  0  in  the  left  plots  of  Fig.  2.  One  can  also  understand  in  this  way  why  the 
difference  disappears  more  slowly  towards  low  E  and  L  along  the  lines 

of  constant  L/E  in  the  upper  left  plot  of  Fig.  2.  The  structures  seen  in  the  plot  emerge 
from  the  matter  dependent  phase  shifts  of  oscillation  probabilities  with  more  or  less  equal 
amplitudes.  This  behavior  can  be  understood  in  the  two- flavor  picture  from  Eq.  (79),  which 
implies  that  p™®-**®''  —  pvacuum  [[nearly  with  E  along  the  lines  of  constant  /S.2v  or, 

equivalently,  constant  LjE.  This  growth  is  modulated  in  the  A2iy-direction  by  the  factor 
sin  A2jy  (sin  A2j/  —  A2y  cos  A21/),  which  describes  nicely  the  details  of  the  structures  extending 
along  lines  of  constant  IS.2y  towards  lower  energies. 

The  right  plots  of  Fig.  2  show  the  probability  difference  —  77““““ |  and  the 

relative  difference  for  the  case  a  =  0.03.  This  hgure  contains  the  structures  stemming  from 
the  atmospheric  resonance  energy,  which  were  already  present  in  the  a  =  0  case  (left  plots). 
In  addition,  similar  structures  show  up  around  the  solar  resonance  energy,  reaching  to  lower 
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Figure  3:  Rows  1,  2,  and  3:  absolute  errors  of  the  three  types  of  expansions  for  Pg/i  as 
functions  of  neutrino  energy  E  and  baseline  length  L,  for  Aml^  =  2  ■  10“^  eV^,  a  =  0.03, 
9i2  =  33°,  023  =  45°,  5cp  =  0,  and  several  values  of  sin^20i3.  Matter  of  constant  density 
p  =  3g/cm^  is  assumed.  Row  4  shows  contours  of  Pg^  calculated  numerically. 


energies  and,  due  to  the  larger  effective  mixing  angle,  also  to  shorter  baselines.  The  right 
plots  of  Fig.  2  show  that  in  the  realistic  case  of  two  mass  squared  differences  matter  effects 
are  rather  important  and  have  to  be  taken  into  account  in  a  large  domain  of  the  physically 
interesting  parameter  space. 

6.2  Comparing  the  accuracy  of  the  three  types  of  expansions 

In  order  to  test  the  accuracy  of  our  analytic  expressions,  we  shall  now  compare  the  values 
of  Pe^  obtained  from  the  expansion  formulas  Eq.  (31)  for  the  double  expansion,  Eqs.  (47) 
and  (48)  for  the  single  expansion  in  a,  and  Eqs.  (65)  and  (66)  for  the  single  expansion  in 
Si3  with  the  exact  values  of  Pg^  calculated  numerically.  Figure  3  shows  the  absolute  errors 

|pexpansion  _  pexact  |  ^  p-g_  relative  errorS  |pexpansion  _  pexact  | ^pexact  three 

types  of  expansions  as  functions  of  the  neutrino  energy  E  and  baseline  length  L  for  various 
values  of  sin^20i3.  For  reference,  we  display  in  the  lowest  row  of  graphs  in  each  of  these 
figures  also  the  probability  Pe^  itself,  calculated  numerically  without  any  approximations. 

First,  we  observe  that  in  general  the  absolute  errors  shown  in  Fig.  3  are  rather  small — 
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Figure  4:  Rows  1,  2,  and  3:  relative  errors  of  the  three  types  of  expansions  for  Pe/x 
fnnctions  of  nentrino  energy  E  and  baseline  length  L,  for  =  2  ■  10“^  eV^,  a  =  0.03, 

6i2  =  33°,  023  =  45°,  5cp  =  0,  and  several  valnes  of  sin^20i3.  Matter  of  constant  density 
p  =  3g/cm^  is  assnmed.  Row  4  shows  contonrs  of  Pg^  calcnlated  nnmerically. 

at  the  0.1%  level  in  a  large  part  of  the  parameter  space,  whereas  the  relative  errors  shown 
in  Fig.  4  are  considerably  larger,  dne  to  the  smallness  of  the  probability  itself.  Next,  we 
note  that  the  series  expansions  in  a,  i.e.,  the  donble  expansion  in  a  and  S13  and  the  single 
expansion  in  a,  are  only  valid  for 

Am^  L  L 

aA  =  —  ^  1 ,  or  — -C  10^  km/GeV ,  (80) 

4P  E 

i.e.,  far  below  the  hrst  solar  maximnm.  The  obvious  reason  is  that  expanding  terms  of  the 
type  sinaA  is  only  valid  if  aA  is  small,  and  hence,  these  two  types  of  expansion  cannot 
account  for  neutrino  oscillations  with  the  “solar”  frequency.  Note  that  for  S13  =  0  the  single 
expansion  in  a  gives  Pg^  =  0,  since  in  that  case  the  lowest-order  term  of  Pg^  is  proportional 
to  a^,  and  our  single  expansion  in  a  only  contains  terms  up  to  hrst  order  in  a.  This  explains 
why  at  very  small  values  of  S13  the  double  expansion  (which  includes  terms  of  second  order 
in  a)  is  more  accurate  than  the  single  expansion  in  a.  In  contrast,  the  neutrino  oscillations 
with  the  “solar”  frequency  (for  which  the  condition  in  Eq.  (80)  is  violated)  are  in  general 
very  well  accounted  for  by  the  single  expansion  in  S13,  since  it  retains  the  exact  dependence 
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of  the  probability  on  a.  This  expansion  is  the  best  one  for  relatively  small  values  of  S13  and 
large  values  of  L/E. 

Figures  3  and  4  allow  us  to  put  the  above  observations  on  a  more  quantitative  basis.  In 
the  region  of  LjE  given  in  Eq.  (80)  (lower-right  parts  of  the  graphs),  where  the  oscillations 
are  mainly  driven  by  the  “atmospheric”  mass  squared  difference  Amg^,  the  double  expansion 
and  the  single  expansion  in  a  work  rather  well.  For  not  too  large  values  of  sin^  26^13  the 
double  expansion  is  most  accurate,  with  absolute  errors  at  the  0.1%  level  and  relative  errors 
not  exceeding  1%  (5%)  for  sin^  26^13  =  0.001  (0.01).  However,  for  values  of  sin^  2^13  close  to 
the  current  upper  bound  the  single  expansion  in  a  becomes  better,  with  relative  errors  smaller 
than  1%.  The  reason  for  this  is  that  for  large  values  of  sin^  2^13  the  neutrino  oscillations  driven 
by  the  “atmospheric”  frequency  and  mixing  angle  613  completely  dominate  the  probability, 
and,  in  addition,  for  relatively  large  values  of  L  and  energies  close  to  10  GeV,  the  atmospheric 
resonance  becomes  important.  Since  the  single  expansion  in  a  is  exact  in  S13,  it  describes  the 
case  of  relatively  large  S13  well,  and  the  atmospheric  resonance  is  also  correctly  accounted 
for.  At  the  same  time,  the  accuracy  of  the  double  expansion  becomes  slightly  worse  near 
the  hrst  atmospheric  maximum  and  the  resonance.  As  was  pointed  out  in  Ref.  [36],  for  large 
values  of  S13  the  accuracy  of  the  double  expansion  can  be  improved  by  replacing  the  term 
proportional  to  in  Eq.  (31)  by  the  term  given  in  Eq.  (47)  as  the  zeroth  order  term  in  the 
single  expansion  in  a.  For  a  =  0  this  term  describes  the  probability  Pe^  exactly  to  all  orders 
in  S13. 

As  can  be  seen  from  Figs.  3  and  4,  the  single  expansion  in  S13  gives  a  rather  poor  de¬ 
scription  of  the  region  of  L/E  dehned  in  Eq.  (80).  The  reason  is  that  with  only  terms  of  hrst 
order  in  S13  it  is  not  possible  to  obtain  a  correct  description  of  neutrino  oscillations  driven 
by  and  ^13.  This  is  also  rehected  by  the  fact  that  the  lowest-order  in  S13  terms  in 

the  eigenvalues  of  the  Hamiltonian  are  C7(s^g)  [see  Eqs.  (27)-(29)].  On  the  other  hand,  the 
single  expansion  in  S13  is  rather  accurate  for  large  values  of  LjE  violating  the  condition  (80) 
(upper-left  parts  of  the  graphs)  and  relatively  small  values  of  ^13:  for  sin^  2^13  <  0.001  the 
accuracy  is  typically  better  than  1%. 

To  conclude  this  subsection,  we  show  in  Fig.  5  which  type  of  expansion  provides  the  most 
accurate  expression  for  Pe^t,  depending  on  the  values  of  the  expansion  parameters  a  and  S13 
and  for  a  number  of  hxed  values  of  neutrino  energy  E  and  baseline  length  L.  These  plots 
change  very  little  when  the  fundamental  neutrino  parameters  Amgj^,  6^12,  6^23,  and  5cp  are 
varied  within  their  allowed  ranges,  and  also  when  one  switches  over  to  the  other  neutrino 
oscillation  channels.  As  expected,  one  observes  as  a  general  trend  that  the  single  expansion 
in  a  is  best  for  small  a  and  large  S13,  whereas  the  single  expansion  in  S13  is  best  for  small 
S13  and  large  a.  The  double  expansion  is  most  accurate  in  a  region  where  a  and  S13  are  of 
comparable  order  of  magnitude. 

In  agreement  with  the  discussion  related  to  Figs.  3  and  4,  we  hnd  from  Fig.  5  that  in  the 
low-energy  regime  E  ~  0.1  GeV  the  single  expansion  in  S13  is  most  accurate  in  almost  the 
entire  a-613  plane.  For  values  of  L/E  satisfying  the  condition  (80)  and  energies  larger  than  a 
few  GeV  the  double  expansion  is  most  accurate  in  a  large  fraction  of  the  a-613  plane;  only  for 
values  of  S13  close  to  the  current  upper  bound  does  the  single  expansion  in  a  become  better. 
Furthermore,  we  note  that  for  E  ~  10  GeV  and  L  >  5500  km  the  atmospheric  resonance  is 
important,  which  leads  to  a  better  accuracy  of  the  single  expansion  in  a  (see  the  rightmost 
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Figure  5;  The  plot  shows  which  type  of  expansion  for  is  most  accnrate,  depending  on 
a  and  S13  and  for  different  values  of  neutrino  energy  E  and  baseline  length  L.  The  values 
of  the  fundamental  neutrino  parameters  used  are  =  2  ■  10“^  eV^,  612  =  33°,  623  =  45°, 

and  (5cp  =  0.  Matter  of  constant  density  p  =  3g/cm^  is  assumed.  The  dotted  lines  indicate 
the  3cr  allowed  range  of  the  parameter  a. 


panel  in  the  third  row  of  Fig.  5). 

We  conclude  that  the  double  expansion  in  both  a  and  S13  is  most  accurate  in  a  wide 
region  of  the  parameter  space,  where  a  and  S13  are  roughly  of  the  same  order  of  magnitude. 
The  single  expansion  in  S13  has  to  be  used  whenever  neutrino  oscillations  with  the  solar 
frequency  are  important  {e.g.,  in  the  low-energy  regime),  whereas  the  single  expansion  in  a 
is  most  accurate  for  values  of  S13  close  to  the  upper  bound,  or  in  cases  where  the  atmospheric 
resonance  is  important. 

6.3  Accuracy  of  the  double  expansion 

Motivated  by  the  fact  that  the  double  expansion  is  most  accurate  in  a  wider  region  of  the 
parameter  space  than  the  single  expansions,  we  present  in  this  subsection  some  more  accuracy 
tests  for  it.  Figures  6  and  7  show  the  relative  errors  for  given  in  Eq.  (31)  and  P^^  given  in 
Eq.  (34),  respectively,  for  neutrinos,  antineutrinos,  and  normal  and  inverted  hierarchies.  We 
note  that  the  relative  accuracy  of  Per  and  1  —  P^e  is  very  similar  to  the  one  of  Pg^  (Fig.  6), 
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Figure  6;  Relative  errors  of  the  double  expansion  for  for  neutrinos  (rows  1  and  2),  and 
P-^  for  antineutrinos  (rows  3  and  4),  for  normal  hierarchy  (NH)  and  inverted  hierarchy  (IH). 
The  values  of  the  fundamental  neutrino  parameters  used  are  lArUg^l  =  2- 10“^  eV^,  a  =  0.03, 
6*12  =  33°,  6*23  =  45°,  (5cp  =  0,  and  several  values  of  sin^  26*13.  Matter  of  constant  density 
p  =  3  g/ cm^  is  assumed. 

and  the  accuracy  of  P^^  and  Prr  is  similar  to  that  of  P^r  (Fig.  7). 

In  the  region  oi  L/E  dehned  in  Eq.  (80),  the  accuracy  of  Pe^  is  roughly  between  1%  and 
5%.  It  should  be  noted  that  the  absolute  errors  of  the  double  expansion  of  Pe^i  depend  very 
weakly  on  the  choice  of  the  fundamental  neutrino  parameters  within  their  allowed  ranges,  and 
they  are  always  very  small,  similar  to  the  hrst  row  of  plots  in  Fig.  3.  However,  the  relative 
errors  shown  in  Fig.  6  are  rather  sensitive  to  the  fundamental  parameters  values  in  the  region 
of  low  L/E.  The  reason  is  that  the  probability  Pe^  itself  is  typically  very  small  in  this  region, 
being  as  tiny  as  10“^,  or  even  smaller.  Hence,  the  relative  error,  being  a  ratio  of  two  small 
numbers,  is  quite  sensitive  to  variations  of  the  parameters  (see,  e.p.,  the  second  column  in 
Fig.  6).  In  contrast,  P^r  is  rather  large,  because  of  the  oscillations  in  “atmospheric”  channel 
driven  by  the  large  and  maximal  or  nearly  maximal  mixing  sin^  26*23  —  1-  This  is  also 

clear  from  Eq.  (34),  since  P^r  has  a  term  of  zeroth  order  in  both  a  and  S13,  while  Pg^  from 
Eq.  (31)  is  of  second  order  in  these  parameters.  Hence,  in  the  case  of  P^ri  the  relative  errors 
are  of  the  same  order  as  the  absolute  errors,  which  are  always  very  small.  The  errors  shown 
in  Fig.  7  are  smaller  than  0.1%  for  sin^  26*13  A  0.01  and  smaller  than  1%  for  sin^  26*13  A  0.1. 
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Figure  7;  Relative  errors  of  the  double  expansion  for  for  neutrinos  (rows  1  and  2),  and 
PjTf  for  antineutrinos  (rows  3  and  4),  for  normal  hierarchy  (NH)  and  inverted  hierarchy  (IH). 
The  values  of  the  fundamental  neutrino  parameters  used  are  lArUg^l  =  2- 10“^  eV^,  a  =  0.03, 
6*12  =  33°,  6*23  =  45°,  (5cp  =  0,  and  several  values  of  sin^  26*13.  Matter  of  constant  density 
p  =  3g/ cm^  is  assumed. 


6.4  Probability  expansions  and  relation  to  experiments 

The  discussion  of  the  accuracy  of  the  different  expansions  can  be  used  to  identify  the  best 
set  of  equations  for  a  given  neutrino  oscillation  experiment.  The  hrst  question  is,  however, 
whether  matter  effects  are  relevant  or  if  the  much  simpler  vacuum  probabilities  can  be  used. 
The  absolute  and  relative  differences  between  matter  and  vacuum  probabilities  vanish,  as 
discussed  in  Sec.  6.1,  for  sufficiently  low  energies  or  for  short  enough  baselines.  The  discussion 
showed  that  in  experimentally  interesting  cases  matter  effects  are  almost  always  relevant  for 
the  ranges  of  L  and  E  that  we  considered.  If  L/P  is  such  that  oscillations  can  be  observed, 
then  for  hxed  LjE  the  matter  effects  are  negligible  only  at  very  small  E  (and  consequently 
small  L).  An  example  where  matter  effects  can  be  ignored  at  the  percent  level  is  given  by 
the  reactor  z/g  disappearance  experiments  with  energies  of  a  few  MeV  and  baselines  up  to  a 
few  kilometers. 

In  order  to  identify  the  best  expansion  for  a  given  experiment,  one  must  take  into  account 
that  the  analyses  of  real  experiments  are  based  on  event  rates  and  not  on  probabilities.  It 
is  useful  to  distinguish  here  between  “point  sources”  (accelerator  beams,  reactors,  the  Sun, 
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supernovae,  etc.)  where  the  observed  neutrinos  come  from  one  point-like  region  in  space,  and 
“extended  sources”  {e.g.,  the  atmosphere)  where  the  observed  neutrinos  come  from  many 
different  positions  in  space.  An  event  rate  based  analysis  of  point  sources  implies  essentially, 
for  a  given  energy,  an  additional  factor  1/L^  in  order  to  account  for  the  reduction  of  the 
flux  as  a  function  of  the  baseline  L.  In  addition,  the  detection  cross  section  is  typically 
proportional  to  E'^,  where  ^  r  <  2,  depending  on  the  energy  range  and  on  the  detection 
process.  In  the  two-flavor  approximation,  when  going  from  the  first  oscillation  maximum  at 
A21/  =  7r/2  to  the  nth  oscillation  maximum  at  A2U  =  (n  —  l/2)7r  this  implies  that  the  event 
rates  drop  by  a  factor  l/(2n  —  1)^  if  the  baseline  is  increased,  or  by  a  factor  l/(2n  —  1)'’  if  the 
energy  is  decreased.  Hence,  e.g.,  for  a  hxed  energy  the  rate  in  the  second  (third)  maximum 
is  already  reduced  by  a  factor  1/9  (1/25).  This  explains  why  most  current  proposals  for  long 
baseline  oscillation  experiments  aim  at  the  first  oscillation  maximum  of  the  “atmospheric” 
oscillations.  The  oscillations  driven  by  the  “solar”  frequency  are  then  a  sub- leading  effect, 
and  in  most  cases  the  double  expansion  works  very  well.  The  overall  relative  precision  of  the 
double  expansion  is  then  typically  in  the  range  of  a  percent,  or  at  worst,  up  to  a  few  percent 
in  this  case.  However,  if  613  is  close  to  its  current  npper  limit,  the  single  expansion  in  a  is 
even  more  precise  than  the  double  expansion.  On  the  other  hand,  if  ^13  is  tiny,  well  below 
the  sensitivity  limits  of  any  planned  accelerator  experiment,  the  single  expansion  in  S13  is 
more  precise  than  the  double  expansion,  especially  for  low  energies.  Note,  however,  that  the 
donble  expansion  has  often  a  snfficient  precision  even  when  it  is  not  the  best  expansion. 

The  discussion  of  point  sonrces  does  not  change  much  for  experiments  which  aim  at 
higher  oscillation  maxima.  The  l/(2n  —  1)^  flux  factor  due  to  the  beam  divergence  and/or 
the  l/(2n  —  1)”  factor  dne  to  the  energy  dependence  of  the  detection  cross  sections  rednce 
the  observed  event  rates  and  limit  all  proposals  to  the  hrst  few  oscillation  maxima  at  best. 
The  double  expansion  still  works  qnite  well  for  snch  proposals,  as  long  as  the  sub-leading 
oscillations  governed  by  the  “solar”  freqnency  stay  in  the  linear  regime,  where  an  expansion 
in  a  makes  sense.  This  leads  to  the  condition  2n  — 1  -C  ~  40,  which  is  fnlhlled  for  the  first 
few  maxima.  Higher  valnes  of  n  are  currently  not  proposed,  since  the  flnx  and/or  cross  section 
wonld  drop  by  a  large  factor,  which  would  make  the  detectors  and  sonrces  nnaffordable.  For 
largest  (smallest)  ^13  the  single  expansion  in  a  (S13)  works  again  nnmerically  even  better  then 
the  donble  expansion,  just  like  in  the  case  of  the  first  oscillation  maximnm. 

The  discussion  is  somewhat  different  for  sonrces  which  are  not  point-like.  In  that  case, 
there  is  in  general  no  1  /L^  suppression  of  the  flux,  and  the  observed  event  rates  are  not  dom¬ 
inated  by  the  first  oscillation  maximum.  An  important  example  are  atmospheric  nentrinos, 
for  which  a  wide  range  of  baselines  and  energies  contribntes:  the  energy  window  of  practical 
interest  ranges  from  a  100  MeV  to  a  few  tens  of  GeV,  while  the  baseline  lies  between  10  km 
and  10^  km.  This  wide  range  of  energies  and  baselines  makes  clear  that  none  of  the  discussed 
expansions  covers  the  fnll  parameter  range  at  the  percent  level.  In  general  the  donble  ex¬ 
pansion  works  qnite  well  for  small  or  moderate  L/E,  while  the  S13  expansion  tends  to  work 
best  for  larger  LjE  values.  At  the  same  time,  close  to  the  atmospheric  resonance  or  if  the 
mixing  angle  ^13  is  near  its  current  npper  bonnd  and  for  not  too  large  L/E,  the  a  expansion 
is  expected  to  be  the  most  relevant  one. 

In  general,  for  all  types  of  experiments  and  for  any  given  valnes  of  L,  E,  and  ^13  the  most 
accnrate  expansion  can  readily  be  fonnd  with  the  help  of  Fig.  5. 
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7  Summary  and  conclusions 


In  this  paper  we  have  presented  three  different  sets  of  approximate  formulas  for  the  prob¬ 
abilities  of  neutrino  oscillations  in  matter  and  in  vacuum.  We  have  shown  that  in  general 
the  probabilities  for  all  possible  oscillation  channels  for  neutrinos  as  well  as  for  antineutrinos 
can  be  obtained  from  just  two  independent  probabilities  by  using  unitarity  and  a  symmetry 
related  to  the  “atmospheric”  mixing  angle  623  [see  Eq.  (14)].  One  possible  choice  for  these 
two  probabilities  is  P^fj,  and  P^r-  We  have  derived  the  expressions  for  the  neutrino  oscilla¬ 
tion  probabilities  in  matter  of  constant  density,  expanded  in  terms  of  the  small  parameters 
a  =  S13  or  in  both  of  them.  Below  we  summarize  the  main  features  of  these 

expansions: 

•  Double  expansion  up  to  second  order  in  ol  and  S13  [Sec.  3.1,  Eqs.  (30)  (34)] 

The  neutrino  oscillation  probabilities  are  expanded  in  both  small  parameters  up  to 
second  order.  In  general,  these  expressions  are  valid  for  oA  =  ^m\^L/{4:E)  1 

OT  L/E  -C  lO^km/GeV,  i.e.,  if  neutrino  oscillations  due  to  the  solar  mass  squared 
difference  Am^i  are  not  important.  The  accuracy  of  the  approximation  is  good  for  a 
wide  range  of  the  parameters.  Typically,  the  relative  errors  of  Per,  and  1  —  Pee  are 
between  1%  and  5%,  the  relative  errors  of  P^r,  and  Prr  and  the  absolute  errors  for 
all  probabilities  are  of  the  order  0.1%. 

•  Single  expansion  up  to  first  order  in  a  [Sec.  4.1,  Eqs.  (45)-(50)] 

The  neutrino  oscillation  probabilities  are  expanded  in  a,  but  the  exact  dependence  on 
Si3  is  retained.  Like  in  the  case  of  the  double  expansion,  these  formulas  are  valid  in 
the  region  where  the  oscillations  driven  by  the  solar  mass  squared  difference  are  not 
important:  a  A  =  Am^iL/ {AE)  or  L/E  lO^km/GeV.  The  accuracy  is  better 
than  the  one  of  the  double  expansion  for  values  of  S13  >  0.1  close  the  current  upper 
bound,  or  if  the  atmospheric  resonance  is  important.  For  instance,  for  a  matter  density 
p  =  3  g/cm^  and  Aml^  =  2  ■  10“^  eV^  this  is  the  case  for  P  ~  10  GeV  and  L  >  5500  km. 

•  Single  expansion  up  to  first  order  in  S13  [Sec.  5.1,  Eqs.  (63)  (68)] 

The  neutrino  oscillation  probabilities  are  expanded  in  S13,  but  the  exact  dependence 
on  a — and  therefore  on  Am^i — is  retained.  Hence,  these  expressions  are  useful  in  the 
region  where  the  oscillations  due  to  the  solar  mass  squared  difference  are  relevant: 
a  A  =  Am2iL/{4E)  >  1,  or  L/E  >  lO^km/GeV.  In  particular,  this  is  the  case  for  low 
energies  E  ~  0.1  GeV  and  L  >  10^  km.  For  values  of  L/E  outside  the  “solar”  regime 
this  type  of  expansion  is  only  useful  for  very  small  values  of  S13  <  10“^  -i- 10“^. 

•  Vacuum  limit  of  the  expansions  [Secs.  3.2,  4.2,  5.2] 

The  vacuum  limit  for  each  type  of  expansion  is  valid  if,  in  addition  to  the  requirements 
ensuring  the  validity  of  a  given  type  of  expansion,  matter  effects  can  be  neglected  (see 
discussion  in  Sec.  6.1).  This  usually  implies  low  energies  or  short  baselines.  We  hnd 
that  in  many  realistic  cases  matter  effects  are  of  the  order  of  a  few  percent.  Hence, 
if  an  accuracy  at  that  level  is  required,  one  should  employ  the  formulas  for  neutrino 
oscillations  in  matter. 
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In  conclusion,  we  have  presented  a  collection  of  formulas  for  three-flavor  neutrino  oscilla¬ 
tion  probabilities  by  deriving  expansions  in  small  parameters.  We  have  performed  a  detailed 
analysis  of  the  accuracy  of  these  expansions  and  determined  the  parameter  regions  where 
they  are  most  accurate.  The  expansions  of  the  neutrino  oscillation  probabilities  in  matter 
of  constant  density  are  useful  for  the  analytical  understanding  of  the  physics  of  future  neu¬ 
trino  oscillation  experiments.  Furthermore,  we  have  also  presented  expansion  formulas  for 
the  neutrino  oscillation  probabilities  in  arbitrary  matter  density  prohles  (see  Appendix  B), 
which  can  be  applied  to  a  large  class  of  problems. 
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A  Hamiltonian  diagonalization  approach 

In  this  appendix,  we  present  the  details  of  the  approximate  diagonalization  of  the  effective 
Hamiltonian  of  the  neutrino  system  in  the  case  of  matter  of  constant  density. 

A.l  Diagonalizing  the  Hamiltonian  perturbatively 

In  order  to  derive  the  double  expansions  given  in  Sec.  3,  we  write  the  Hamiltonian  of  Eq.  (6) 
as 

H^^OMMUlOg,  (Al) 

where  M  =  O13O12  diag  (0,  a,  1)  0^0^  -|-  diag(A,0,0),  and  the  matrices  Oij  and  Us  have 
been  dehned  after  Eq.  (1).  The  matrix  M  can  be  explicitly  written  as  M  =  {Amh/2E)~^H' 
by  setting  (5cp  =  0  in  77',  which  is  given  in  Eq.  (Bl)  below.  First,  we  diagonalize  the  matrix 
M  hj  M  =  WMW^  with  M  =  diag(Ai,  A2,  A3)  and  W  being  a  unitary  diagonalizing  matrix. 
This  diagonalization  is  performed  by  using  perturbation  theory  up  to  second  order  in  the 
small  parameters  a  and  S13,  he.,  we  write  M  =  where 

contains  all  terms  of  hrst  (second)  order  in  a  and  S13.  One  hnds 

=  diag(A,  0, 1)  =  diag(AS°\  A^°\  A^°^) ,  (A2) 
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M(2) 


as?2 

OiSi2Ci2  Si3  \ 

OiSi2Ci2 

ac?2  0 

Sl3 

0  0  / 

*13 

0  -asi3S?2 

0 

0  — Q;<S13Si2Ci2 

~Q:*13*i2 

— Q;Si3Si2Ci2  ~*13 

(A3) 

(A4) 


For  the  eigenvectors  we  write  n*  =  n®  since  is  diagonal  at  zeroth  order, 

we  have  vf'^  =  Cj.  Then,  the  first  and  second  order  corrections  to  the  eigenvalues  are  given 
by 


a:'> 

Af> 


M, 


(1) 


(2) 


E 


x(o)  x{o)  ’ 

jV*  A' 


and  the  corrections  to  the  eigenvectors  are  calculated  by 


= 


(2) 

v]  = 


E 


(1) 


y(o)  _  ^(0)  ’ 


j^i  ''j 

1 


E 


xW  x(o) 

jj^i  \  A 


Mjf  +  _  a(«  („w) 


Cj. 


(AS) 

(A6) 

(A7) 

(AS) 


The  mixing  matrix  in  matter  is  given  by  17'  =  O23USW  with  W  =  {vi,V2,V3),  and  the 
eigenvalues  of  the  Hamiltonian  [see  Eqs.  (27)-(29)]  are  obtained  as  Ei  =  [Amh/ {2E)]\i 
(*  =  1,2,3). 

In  our  paper,  we  do  not  in  general  order  the  eigenvalues  according  to  their  magnitude. 
Such  an  ordering  would  create  problems  as  one  would  have  to  re-label  the  eigenvalues  upon 
passing  through  each  of  the  two  MSW  resonances.  The  ordering  is  actually  unimportant  if 
one  is  careful  to  assign  the  correct  eigenvector  to  each  eigenvalue. 

We  also  reiterate  the  point  discussed  at  the  end  of  Sec.  3.1:  despite  the  fact  that  in  the 
case  of  the  double  expansion  the  eigenvalues  as  well  as  certain  entries  of  the  leptonic  mixing 
matrix  in  matter  are  divergent  at  A  —a  0  and  A  —a  1,  the  neutrino  oscillation  probabilities 
are  hnite  in  these  limits.  In  particular,  the  correct  vacuum  probabilities  are  recovered  in  the 
limit  A  —A  0.  The  mentioned  divergences  are  of  no  concern  to  us,  since  we  are  interested 
in  oscillation  probabilities  rather  than  in  eigenvalues  or  the  mapping  between  the  mixing  in 
matter  and  in  vacuum,  as  was  the  case,  e.g.,  in  Ref.  [36]. 


A. 2  The  Cayley-Hamilton  formalism 


In  order  to  derive  a  series  expansion  for  the  evolution  matrix  S  in  small  parameters,  we  can 
use  the  following  formula  from  Refs.  [12, 14]: 

S{t,to)  = 

g-iElL 

=  - - - - -  \E‘2E^\2,  ~  {^2  +  E‘^)H  +  H‘^'\ 

{E^  -  E2){Ei  -  E^)  L  2  3  3  I  2W  3;  w  J 

ii?2i 

+  7 - 77 - 7  [£^1^313  -  {El  +  E^)H  + 

{E2-Ei){E2-E^)^  I  3;  w  J 

iBsL 

+  77^ - ^7-77^ - —\EiE2l^-{Ei  +  E2)H  +  H^],  (A9) 

{E^  -  Ei){E:i  -  E2)  ^  1  2  3  V  1  2>  J,  V  ; 

where  H  is  the  Hamiltonian  in  the  flavor  basis  as  given  in  Eq.  (6)  and  Ei  {i  =  1,2,3)  are 
the  eigenvalues  of  this  Hamiltonian.  Inserting  series  expansions  of  the  eigenvalues  [see,  e.g.^ 
Eqs.  (40)-(42)  and  (58)-(60)]  and  the  Hamiltonian  H  into  Eq.  (A9),  it  is  straightforward 
to  obtain  an  expression  for  Sag-,  and  subsequently,  for  the  neutrino  oscillation  probabilities 
Pag  =  \Sg^\\ 


B  Perturbative  expansion  of  the  evolution  equation  for  arbitrary 
matter  density  profiles 


In  this  appendix,  we  present  the  details  of  the  perturbative  expansion  of  the  neutrino  evolution 
equation  in  the  case  of  matter  with  an  arbitrary  density  profile.  We  also  give  the  formulas 
relevant  for  the  particular  case  of  matter  of  constant  density. 

Consider  the  neutrino  evolution  matrix  S'(t,  to)  defined  in  Eq.  (9).  It  is  convenient  to 
rotate  it  by  the  angle  ^23  according  to  Eq.  (13).  The  evolution  matrix  in  the  rotated  ba¬ 
sis  S"(t,  to)  satisfies  the  evolution  equation  i(d/dt)5"(t,  to)  =  H'{t)S'{t,to)  with  the  initial 
condition  S' {to, to)  =  1  and  the  Hamiltonian  H'{t)  given  by  [28] 


H'{t) 


Here 


'512^13^21  +  -Sio A31  -|-  V (t) 

’S12C12C13A21 
SisCis  (A31  —  5^2^21) 


’S12C12C13A21  S13C13  (A31  —  S^2^2i)  e 

C12A21  — "Si2Ci2Si3e  ‘'^'^^A21 

~'Si2Ci2Si3e^^‘^^A2i  C^3A31  +  S^2'5l3A21 


2E 


A31 


2E 


It  is  important  to  notice  that  H'{t)  [and  hence,  S'{t,  to)]  does  not  depend  on  the  mixing  angle 
623.  This  is  a  consequence  of  the  specific  parameterization  of  the  leptonic  mixing  matrix  U 
in  Eq.  (1),  in  which  the  matrix  O23  is  the  leftmost  one,  and  the  fact  that  the  matrix  of  the 
matter-induced  potentials  in  the  effective  Hamiltonian  H{t)  commutes  with  O23.  Once  the 
matrix  S'{t,to)  is  found,  one  can  use  Eq.  (13)  to  rotate  back  to  the  original  flavor  basis. 
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It  is  convenient  to  decompose  the  Hamiltonian  H'{t)  as  H'{t)  =  H^lt)  +  Hj,  where  -Ho(^) 
is  of  zeroth  order  in  some  expansion  parameter  and  H'j  is  the  remainder.  Accordingly,  the 
evolution  matrix  can  be  written  as  S'(t,to)  =  SQ(t,to)S[(t,to),  where  So(t,to)  satishes 

fo)  =  i/o(t)Ao(t,  to) ,  *S'o(to,  to)  =  1 .  (B3) 

Then  the  matrix  S'[(t,to)  satisfies 

i^^i(t,to)  =  [A'(t,to)-'hr;^^(t,to)]^U^,^o),  ^U^o,to)  =  1.  (B4) 

Up  to  now,  no  approximations  have  been  made.  Next,  we  shall  find  the  evolution  matrix 
perturbatively.  Let  H[  be  the  part  of  Hj  that  is  of  the  first  order  in  the  chosen  expansion 
parameter.  Then,  to  first  order  in  this  parameter,  one  hnds 

S'{t,  to)  ~  ^^(t,  to)  -  i^^(t,  to)  Aa' (t',  to)]  dt' ,  (B5) 

Jto 

from  which,  upon  rotation  back  to  the  original  flavor  basis,  the  evolution  matrix  5'(t,to)  is 
obtained.  We  shall  now  consider  the  application  of  the  above  formalism  to  the  cases  when 
the  expansion  parameter  is  either  a  or  S13. 

B.l  Expansion  up  to  first  order  in  cx 

In  this  case,  the  zeroth-order  Hamiltonian  iLo(^)  obtained  from  Eq.  (Bl)  by  taking  the 
limit  a  — 0  (he.,  A21  — 0).  The  zeroth-order  evolution  matrix  SQ{t,to)  can  then  be  written 
as  [26] 

/  u(t,to)  0  \ 

^'(t,to)=  0  fit, to)  0  ,  (B6) 

\  -v*it,to)e^^^^  0  u*it,to)  ) 

where  uit,  to)  and  vit,  to)  are  to  be  found  from  the  solution  to  the  two-flavor  problem  governed 
by  the  mass  squared  difference  Amh,  mixing  angle  6*13,  and  matter  potential  Vit),  and 

fit,  to)  =  exp  ^i  Ait')  df'J  ,  Ait)  ^  +  U(f)J  .  (B7) 

The  parameters  u  and  v  satisfy  |mP  +  |np  =  1.  The  factor  was  pulled  out  from  v  for 

convenience  (this  will  simplify  some  of  the  formulas  given  below). 

The  neutrino  oscillation  probabilities  to  hrst  order  in  a  in  the  case  of  matter  with  an 
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arbitrary  density  profile  can  be  written  as 


-Pp.fi  =  \u\ 


2  — —  (2  cos  26*13  -^1  +  sin  26*13  -^2)] , 

Jj 

aA 


(B8) 


aA 


Pefj,  =  'S23|'yr  +  2  — ^  -5^2  Ini[nn*  (2  cos  26*13  -^1  +  sin  26*13  -^2)]  H — ^  sin  26*12  sin  26*23 

1^  Jj 

X  {cos  (5cp  Ini[n*/*(ci3-^3  +  -Si3-^4)]  +  sin  (5cp  Re[n*/*(ci3-^3  +  -S13/4)]}  ,  (B9) 


P^r  =  7  sin^  26*23  [1  +  “  2Re(M/)] 


aA 


sin^  26*23  s  s^2 


+ 


(m  -  f*)v*  ( cos  26*13  -^1  +  ^  sin  26*13  h 


C12P  —  S^2  (  '^13  -^5  +  Ci3  Iq  —  -  sin  2^13  Jy 


Im(n/) 


aA 


-  —  sin  26*12  sin26*23|  cos(5cp[c23lni[(n  - /*)n*(ci3/3  +  S13/4)] 

+  Ini[(c23|Mp  -  S23  -  cos  26*23  m/)  (C13/9  -  sis/s)]] 

+  sin  (5cp  [c23  Re[(M  -  /*)n*(ci3  h  +  Si3  -^4)] 

-  Re[(c23|n|^  -  S23  -  cos 26*23  uf)  (C13/9  -  sis/g 
where  the  qnantities  Ii,  I2, . . . ,  I9  are  given  by 

h=  [  u*{t',to)v{t',to)dt' ,  l2=  [  [v{t',tof -u*{t',tof]dt' . 


'to 


'to 


(BIO) 


h=  u*{t',to)f{t',to)dt\  h=  v{t',to)f{t',to)dt', 


'to 


'to 


h=  \u{t',to)\‘^  dt' , 


h  =  i  \v{t\tQ)\^  dt\ 


'to 


'to 


Jy  =  /  [u{t\to)v{t\tQ)  +U*(t\tQ)V*{t\tQ)]dt' , 


'to 


h=  u*{t',to)f*{t',to)dt' ,  h  =  i  v{t\tQ)f*{t\tQ)dt' . 


'to 


'to 


Note  that  the  integrals  Is,  Jg,  and  are  real,  Jg  +  Jg  =  t  —  tg,  and  the  integrals  Jg  and  Jg  can 
be  obtained  from  and  1^,  respectively,  by  snbstitnting  /  — f*. 

In  the  case  of  matter  of  constant  density,  the  parameters  u,  v,  and  /  are  given  by 


Ci^At  .cos 2^13  — A  .  Ci^At 
U[t,  0)  =  cos - - - h  1 - y:: - Sm 


C 


13 


.sin 26*13  .  Ci^At 

V{t,  0)  =  — 1 - yy -  sm  ■ 


C 


(Bll) 

(B12) 
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f{t,  0)  =  exp 


i(i+^i 


(B13) 


Inserting  these  expressions  into  the  integrals  Ii,  I2, . . . ,  I9  and  the  obtained  results  into  Eqs.  (B8)- 
(BIO)  leads  to  Eqs.  (45)-(50). 


B.2  Expansion  np  to  first  order  in  S13 


This  case  was  studied  in  detail  in  Ref.  [28]  (Appendix  A),  which  we  closely  follow  here.  The 
zeroth-order  Hamiltonian  is  obtained  from  Eq.  (Bl)  by  taking  the  limit  ^13  — 0.  The 

zeroth-order  evolution  matrix  S'Q(f,to)  can  be  written  as 

/  x{t,to)  y{t,tQ)  0  \ 

5'o(^,fo)=  x*{t,to)  0  ,  (B14) 

V  0  0  9{t,to)  J 


where  x(t,  to)  and  y(t,  to)  are  to  be  found  from  the  solution  to  the  two-flavor  problem  governed 
by  the  mass  squared  difference  Amh,  mixing  angle  6^12,  and  matter  potential  V(t),  and 


g{t,to)  =  exp 


— i  f  A{t')  dt'  , 


A(t) 


2E 


1 

2 


Aml^ 

2E 


+  V{t) 


(B15) 


The  parameters  x  and  y  satisfy  |a;p-|-  \y\‘^  =  1. 

The  neutrino  oscillation  probabilities  to  first  order  in  S13  in  the  case  of  matter  with  an 
arbitrary  density  prohle  can  be  written  as 

Pee=\x\^,  (B16) 

Pe^,  =  C23  \y?  -  sin  26123  ,  (B17) 

=  S23C23  k  -  9*?  +  sin  2^23  Im[(a;  -  9*)9{sloB  -  cl^D)] ,  (B18) 


where  the  quantities  i?,C,  and  D  are  dehned  as  [47] 


B  =  B(t,  to)  =  aly*^t(t,  to)  +  blx,t(t,  to) ,  (B19) 

C  ^  C{t,to)  =  -  b*i;,^{t,to) ,  (B20) 

D  =  D{t,to)  =  a*I**^t^{t,to)  +  b*Ilt^{t,to) .  (B21) 


Here  the  integral  is  dehned  as 


I^,sit,to)=  (p{t',s)9{t\s)dt' , 
Jto 


and  the  quantities  a  and  b  are  given  by 

Am? 


31 


2E 


b  = 


2E 


(B22) 


■Sl3(l-4«)e-'^°^ 

(B23) 

2 

21  — i(5pp 

— -Sl3Si2Ci2e  CP  . 

(B24) 
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In  the  case  of  matter  of  constant  density,  the  parameters  x,  y,  and  g  are  given  by 


aCu^t  cos  26^12  —  A/o;  Q;Ci2At 
x{t,  0)  =  cos - - - h  1 - - sm  ■ 


C 


12 


.sin 26*12  .  aCuAt 
=  — 1 - sm  ■ 


g{t,  0)  =  exp 


Ci2  L 

.{A  +  a-2)A 


(B25) 

(B26) 

(B27) 


Inserting  these  expressions  into  Eqs.  (B19)  -  (B22)  and  the  obtained  results  into  Eqs.  (B16)- 
(B18)  leads  to  Eqs.  (63)-(68). 
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